Problem description¶
In the telecom industry, customers are able to choose from multiple service providers and actively switch from one operator to another. In this highly competitive market, the telecommunications industry experiences an average of 15-25% annual churn rate. Given the fact that it costs 5-10 times more to acquire a new customer than to retain an existing one, customer retention has now become even more important than customer acquisition.
For many incumbent operators, retaining high profitable customers is the number one business goal. To reduce customer churn, telecom companies need to predict which customers are at high risk of churn.
Definition of churn¶
There are two main models of payment in the telecom industry - postpaid (customers pay a monthly/annual bill after using the services) and prepaid (customers pay/recharge with a certain amount in advance and then use the services).
In the Indian and the southeast Asian market, approximately 80% of revenue comes from the top 20% customers (called high-value customers). Thus, if we can reduce churn of the high-value customers, we will be able to reduce significant revenue leakage.
Goal of the analysis¶
Define high-value customers based on a certain metric (mentioned later below) and predict churn only on high-value customers.
Import libraries and load datafile¶
In [1]:
# Import required libraries
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
# Until fuction: line seperator
def print_ln():
print('-'*80, '\n')
pd.options.display.float_format = '{:.2f}'.format
# Load csv data file
telecom_orig = pd.read_csv('../data/raw/telecom_churn_data.csv', low_memory=False)
Checking the shape¶
In [2]:
# Understand the properties of loaded dataframe
telecom = telecom_orig.copy()
print('Dataframe Shape: ', telecom.shape); print_ln();
print("Dataframe Info: \n"); telecom.info(); print_ln();
telecom.head(5)
Out[2]:
Checking Categorical and Numerical Columns¶
In [3]:
def type_features(data):
categorical_features = data.select_dtypes(include = ["object"]).columns
numerical_features = data.select_dtypes(exclude = ["object"]).columns
print( "categorical_features :",categorical_features)
print_ln()
print("numerical_features:",numerical_features)
print_ln()
return categorical_features, numerical_features
telecom_cat_features, telecom_num_features = type_features(telecom)
Missing Vlaue Analysis¶
- Delete: Delete the missing values
- Impute:
- Imputing by a simple statistic: Replace the missing values by another value, commonly the mean, median, mode etc.
- Predictive techniques: Use statistical models such as k-NN, SVM etc. to predict and impute missing values
In [4]:
# summing up the missing values (column-wise)
telecom.isnull().sum()
Out[4]:
Percentage of missing values (column-wise)¶
In [5]:
round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
Out[5]:
In [6]:
# Columns with more than 60% missing values
colmns_missing_data = round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
colmns_missing_data[colmns_missing_data >= 60]
Out[6]:
Treating Categorical Variables¶
There are 8 of them and all of these have more than 60% missing values. Missing values can be imputed with '-1' where -1 can be treated as a new category.
CATEGORICAL VARIABLES:
- night_pck_user_6
- night_pck_user_7
- night_pck_user_8
- night_pck_user_9
- fb_user_6
- fb_user_7
- fb_user_8
- fb_user_9
In [7]:
# Missing values can be imputed with '-1'
# where -1 can be treated as a new category.
categorical_vars = ["night_pck_user_6", "night_pck_user_7", "night_pck_user_8",
"night_pck_user_9", "fb_user_6", "fb_user_7", "fb_user_8", "fb_user_9"]
for var in categorical_vars:
telecom[var].fillna(value=-1, inplace=True)
In [8]:
# Create dummy variables for all the six categorial variables
night_pck_user_6_s = pd.get_dummies(telecom['night_pck_user_6'],prefix="night6",drop_first=True)
night_pck_user_7_s = pd.get_dummies(telecom['night_pck_user_7'],prefix="night7",drop_first=True)
night_pck_user_8_s = pd.get_dummies(telecom['night_pck_user_8'],prefix="night8",drop_first=True)
fb_user_6_s = pd.get_dummies(telecom['fb_user_6'],prefix="fb6", drop_first=True)
fb_user_7_s = pd.get_dummies(telecom['fb_user_7'],prefix="fb7",drop_first=True)
fb_user_8_s = pd.get_dummies(telecom['fb_user_8'],prefix="fb8",drop_first=True)
telecom = pd.concat([telecom,night_pck_user_6_s],axis=1)
telecom = pd.concat([telecom,night_pck_user_7_s],axis=1)
telecom = pd.concat([telecom,night_pck_user_8_s],axis=1)
telecom = pd.concat([telecom,fb_user_6_s],axis=1)
telecom = pd.concat([telecom,fb_user_7_s],axis=1)
telecom = pd.concat([telecom,fb_user_8_s],axis=1)
# drop all categorial columns
telecom = telecom.drop(['fb_user_6', 'fb_user_7', 'fb_user_8', 'night_pck_user_6', 'night_pck_user_7', 'night_pck_user_8'], 1)
In [9]:
# Columns with more than 70% missing values
colmns_missing_data = round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
colmns_missing_data[colmns_missing_data >= 70]
Out[9]:
Data Preparation¶
Below Columns used for filtering high-value customers, can't be dropped
- total_rech_data_6
- total_rech_data_7
- av_rech_amt_data_6
- av_rech_amt_data_7
All these columns has more than 70% missing values
In [10]:
print(telecom.total_rech_data_6.describe()); print_ln()
print(telecom.total_rech_data_7.describe()); print_ln()
print(telecom.av_rech_amt_data_6.describe()); print_ln()
print(telecom.av_rech_amt_data_7.describe()); print_ln()
In [11]:
# Let us look at the distribution to find out the value to fill
# total_rech_data_6
plt.figure(1)
fig, ax = plt.subplots(figsize=(10,10))
plt.subplot(2,2,1)
plt.title("total_rech_data_6")
plt.hist(telecom.total_rech_data_6.dropna(),bins=25)
# total_rech_data_7
plt.subplot(2,2,2)
plt.title("total_rech_data_7")
plt.hist(telecom.total_rech_data_7.dropna(),bins=25)
# av_rech_amt_data_6
plt.subplot(2,2,3)
plt.title("av_rech_amt_data_6")
plt.hist(telecom.av_rech_amt_data_6.dropna(),bins=25)
# total_rech_data_7
plt.subplot(2,2,4)
plt.title("av_rech_amt_data_7")
plt.hist(telecom.av_rech_amt_data_7.dropna(),bins=25)
Out[11]:
Insights from the distribution graphs¶
- Variables are numeric.
- The minimum value is 1, and not zero,
In [12]:
# missing -- indicating that no recharge was done in this month.
# Thus, missing values can be imputed with 0.
telecom.total_rech_data_6.fillna(value=0, inplace=True)
telecom.total_rech_data_7.fillna(value=0, inplace=True)
telecom.total_rech_data_8.fillna(value=0, inplace=True)
telecom.total_rech_data_9.fillna(value=0, inplace=True)
telecom.av_rech_amt_data_6.fillna(value=0, inplace=True)
telecom.av_rech_amt_data_7.fillna(value=0, inplace=True)
telecom.av_rech_amt_data_8.fillna(value=0, inplace=True)
telecom.av_rech_amt_data_9.fillna(value=0, inplace=True)
In [13]:
#Column analysis % of missing values
print(round(100*(telecom.isnull().sum()/len(telecom.index)), 2)); print_ln()
# Columns with more than 70% missing values
colmns_missing_data = round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
print(colmns_missing_data[colmns_missing_data >= 70])
In [14]:
drop_missing_values_columns = colmns_missing_data[colmns_missing_data > 70].index
telecom.drop(drop_missing_values_columns, axis=1, inplace=True)
print('Number of columns droped here = ', len(drop_missing_values_columns))
print('Dataframe Shape: ', telecom.shape)
Filter High-Value Customers¶
In [15]:
# New DERIVED Columns for Total Data Recharge Amount
telecom["total_rech_data_amt_6"] = telecom["total_rech_data_6"]*telecom['av_rech_amt_data_6']
telecom["total_rech_data_amt_7"] = telecom["total_rech_data_7"]*telecom['av_rech_amt_data_7']
telecom["total_rech_data_amt_8"] = telecom["total_rech_data_8"]*telecom['av_rech_amt_data_8']
telecom["total_rech_data_amt_9"] = telecom["total_rech_data_9"]*telecom['av_rech_amt_data_9']
# We can drop total_rech_data_* and av_rech_amt_data_*
drop_data_columns = ["total_rech_data_6", "total_rech_data_7", "total_rech_data_8", "total_rech_data_9",
'av_rech_amt_data_6', 'av_rech_amt_data_7', 'av_rech_amt_data_8', 'av_rech_amt_data_9']
telecom.drop(drop_data_columns, axis=1, inplace=True)
In [16]:
# Filter high-value customers
avg_rech_amt_6_7 = (telecom["total_rech_amt_6"] + telecom["total_rech_data_amt_6"] +
telecom["total_rech_amt_7"] + telecom["total_rech_data_amt_7"]) / 2.0
amont_70_pc = np.percentile(avg_rech_amt_6_7, 70.0)
print('70 percentile of first two months avg recharge amount: ', amont_70_pc)
telecom = telecom[avg_rech_amt_6_7 >= amont_70_pc]
print('Dataframe Shape: ', telecom.shape)
Identify Curn & Tag Churners¶
In [17]:
# Identify Churn
X = telecom["total_ic_mou_9"] + telecom["total_og_mou_9"] + telecom["vol_2g_mb_9"] + telecom["vol_3g_mb_9"]
telecom["churn"] = np.where(X, 0, 1)
#telecom["churn"].head(10)
# Columns to be dropped: all columns ending with "_9"
drop_fourth_month_columns = [hdr for hdr in list(telecom) if hdr.endswith("_9")]
drop_fourth_month_columns += ["sep_vbc_3g"]
telecom.drop(drop_fourth_month_columns, axis=1, inplace=True)
print('Number of columns droped here = ', len(drop_fourth_month_columns))
print('Dataframe Shape: ', telecom.shape); print_ln();
Identify columns that have no varience & Drop¶
In [18]:
#Identify columns that have no varience
telecom_unique_count = telecom.nunique().sort_values(ascending=False)
# Identify bad colums that has no infromation (all entries are NA or same)
# Find columns with all NULL entries and add to drop_columns list
telecom_unique_count_is_zero = telecom_unique_count[telecom_unique_count == 0]
print("Dataframe Unique Value Count is ZERO (all null values): \n", telecom_unique_count_is_zero); print_ln();
drop_no_varience_columns = list(telecom_unique_count_is_zero.index)
# Find columns with all same entries and add to drop_columns list
telecom_unique_count_is_one = telecom_unique_count[telecom_unique_count == 1]
print("Dataframe Unique Value Count is ONE (all same values): \n", telecom_unique_count_is_one); print_ln();
drop_no_varience_columns += list(telecom_unique_count_is_one.index)
telecom.drop(drop_no_varience_columns, axis=1, inplace=True)
print('Dataframe Shape: ', telecom.shape); print_ln();
print('Number of columns dropped here = ', len(drop_no_varience_columns))
print('Dataframe Shape: ', telecom.shape); print_ln();
ROW-WISE: MISSING VALUES¶
In [19]:
# sum it up to check how many rows have all missing values
print("Rows with all NULL values =", telecom.isnull().all(axis=1).sum())
# sum of misisng values in each row
rows_missing_data = telecom.isnull().sum(axis=1)
rows_missing_data[rows_missing_data > 0]
Out[19]:
In [20]:
#small percentage of rows has more than or equal to 27 missing values
drop_rows = telecom.isnull().sum(axis=1) >= 27
print("Number of rows with >= 27 missing values = {}".format(drop_rows.sum()))
print(" % of rows with >= 27 missing values = {} %".format(drop_rows.sum()/telecom.shape[0]*100))
# drop rows with missing values greater >= 27
telecom = telecom[telecom.isnull().sum(axis=1) < 27]
print('Dataframe Shape after dropping rows: ', telecom.shape); print_ln();
In [21]:
rows_missing_data = telecom.isnull().sum(axis=1)
print(rows_missing_data[rows_missing_data > 0])
In [22]:
# look at the COLUMN summary again
df = round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
df[df > 0]
Out[22]:
The remaining missing values are in date columns, before imputing convert them into appropreate values for analysis.
Fix Dates & convert to meaningful numbers¶
In [23]:
# covert date columns to python datetime format
date_vars = ["date_of_last_rech_6", "date_of_last_rech_7", "date_of_last_rech_8"]
telecom[date_vars].head()
Out[23]:
In [24]:
# convert to datetime
for col in date_vars:
telecom[col] = pd.to_datetime(telecom[col])
print(telecom[date_vars].info())
telecom[date_vars].head()
Out[24]:
In [25]:
# Create new days columns, instead of date
import datetime
last_date_of_month_6 = datetime.datetime.strptime("30-06-2014", "%d-%m-%Y")
last_date_of_month_7 = datetime.datetime.strptime("31-07-2014", "%d-%m-%Y")
last_date_of_month_8 = datetime.datetime.strptime("31-08-2014", "%d-%m-%Y")
telecom["rech_days_left_6"] = (last_date_of_month_6 - telecom.date_of_last_rech_6).astype('timedelta64[D]')
telecom["rech_days_left_7"] = (last_date_of_month_7 - telecom.date_of_last_rech_7).astype('timedelta64[D]')
telecom["rech_days_left_8"] = (last_date_of_month_8 - telecom.date_of_last_rech_8).astype('timedelta64[D]')
day_columns = ["rech_days_left_6", "rech_days_left_7", "rech_days_left_8"]
print(telecom[day_columns].info())
In [26]:
# Drop all old date columns: add dates columns to drop_column list
telecom.drop(date_vars, axis=1, inplace=True)
telecom1 = telecom.copy()
print('Number of columns droped here = ', len(date_vars))
print('Dataframe Shape: ', telecom.shape)
In [27]:
# look at the COLUMN summary again
df = round(100*(telecom.isnull().sum()/len(telecom.index)), 2)
df[df > 0]
Out[27]:
In [28]:
telecom[['rech_days_left_6', 'rech_days_left_7', 'rech_days_left_8']].describe()
Out[28]:
In [29]:
# Mean & Median are different: imputing with median value for recharge days left
telecom.loc[np.isnan(telecom['rech_days_left_6']), ['rech_days_left_6']] = telecom['rech_days_left_6'].median()
telecom.loc[np.isnan(telecom['rech_days_left_7']), ['rech_days_left_7']] = telecom['rech_days_left_7'].median()
telecom.loc[np.isnan(telecom['rech_days_left_8']), ['rech_days_left_8']] = telecom['rech_days_left_8'].median()
print("\nTotal Number of missing values in the data frame = {}".format(telecom.isnull().sum().sum()))
print('Dataframe Shape: ', telecom.shape); print_ln();
Outlier Assessment¶
Use data distribution to find outliers
In [30]:
# Checking outliers at 25%,50%,75%,90%,95% and 99%
telecom.describe(percentiles=[.25,.5,.75,.90,.95,.99])
Out[30]:
- Most of the outliers are on the revenue data and the model being built is on usage of voice and data.Hence there is no need to treat these outliers
In [31]:
print('Original Dataframe Shape: ', telecom_orig.shape); print_ln();
print("Original Dataframe Info: \n"); telecom_orig.info(); print_ln();
print("Original Dataframe Nulls:", telecom_orig.isnull().sum().sum()); print_ln();
print('Cleaned Dataframe Shape: ', telecom.shape); print_ln();
print("Cleaned Dataframe Info: \n"); telecom.info(); print_ln();
print("Cleaned Dataframe Nulls:", telecom.isnull().sum().sum()); print_ln();
Save the cleaned data in new file¶
In [32]:
# write treated telecom file
telecom.to_csv("../data/processed/telecom_churn_data_clean.csv", sep=',', index=False)
Dependent Variable: check data imbalance¶
In [33]:
#checking whether data is imbalanced or not
churn = telecom["churn"].value_counts()
churn_rate = churn[1]/(churn[1] + churn[0])
print("Total Churn Count = {}".format(churn[1]))
print("Total Not-Churn Count = {}".format(churn[0]))
print("Churn Rate = {:.2f}%".format(churn_rate*100))
print_ln()
Data is highly imbalanced, we need to balance it before we can create the model
Data Analysis¶
1) Univariate Analysis¶
- now use histogram or probability density function for univiraite analysis
- one variable analysis
- this is important to know feature importance in classification
1.1) Univariate Analysis on 'total_rech_num_6', 'total_rech_num_7', 'total_rech_num_8¶
In [34]:
#Univariate Analysis on 'total_rech_num_6', 'total_rech_num_7', 'total_rech_num_8
df_not_churn = telecom[telecom["churn"] == 0]
df_churn = telecom[telecom["churn"] == 1]
In [35]:
#claculation of mean for churn or non churn
print("Mean:")
print(np.mean(df_not_churn["total_rech_num_6"]))
print(np.mean(df_churn["total_rech_num_6"]))
print("median:")
print(np.median(df_not_churn["total_rech_num_6"]))
print(np.median(df_churn["total_rech_num_6"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_num_6"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_num_6"],np.arange(0,100,25)))
In [36]:
#feature importance of total_rech_num_6
sns.set_style('whitegrid')
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_num_6").add_legend()
#its difficult to know about the churning rate based on total_rech_num_6 feature as overlaped data.
Out[36]:
In [37]:
print("Mean:")
print(np.mean(df_not_churn["total_rech_num_7"]))
print(np.mean(df_churn["total_rech_num_7"]))
print("median:")
print(np.median(df_not_churn["total_rech_num_7"]))
print(np.median(df_churn["total_rech_num_7"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_num_7"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_num_7"],np.arange(0,100,25)))
In [38]:
#feature importance of total_rech_data_7
sns.set_style('whitegrid')
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_num_7").add_legend()
#its difficult to know about the churning rate based on total_rech_data_7 feature as overlaped data.
Out[38]:
In [39]:
print("Mean:")
print(np.mean(df_not_churn["total_rech_num_8"]))
print(np.mean(df_churn["total_rech_num_8"]))
print("median:")
print(np.median(df_not_churn["total_rech_num_8"]))
print(np.median(df_churn["total_rech_num_8"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_num_8"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_num_8"],np.arange(0,100,25)))
In [40]:
#feature importance of total_rech_amt_8
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_num_8").add_legend()
Out[40]:
In [41]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_num_6","total_rech_num_7").add_legend()
Out[41]:
In [42]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_num_7","total_rech_num_8").add_legend()
Out[42]:
In [43]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_num_6","total_rech_num_8").add_legend()
Out[43]:
1.2) Univariate Analysis on 'total_rech_amt_6', 'total_rech_amt_7', 'total_rech_amt_8'¶
In [44]:
print("Mean:")
print(np.mean(df_not_churn["total_rech_amt_6"]))
print(np.mean(df_churn["total_rech_amt_6"]))
print("median:")
print(np.median(df_not_churn["total_rech_amt_6"]))
print(np.median(df_churn["total_rech_amt_6"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_amt_6"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_amt_6"],np.arange(0,100,25)))
In [45]:
#feature importance of total_rech_num_6
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_amt_6").add_legend()
Out[45]:
In [46]:
print("Mean:")
print(np.mean(df_not_churn["total_rech_amt_7"]))
print(np.mean(df_churn["total_rech_amt_7"]))
print("median:")
print(np.median(df_not_churn["total_rech_amt_7"]))
print(np.median(df_churn["total_rech_amt_7"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_amt_7"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_amt_7"],np.arange(0,100,25)))
In [47]:
#feature importance of total_rech_num_6
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_amt_7").add_legend()
Out[47]:
In [48]:
print("Mean:")
print(np.mean(df_not_churn["total_rech_amt_8"]))
print(np.mean(df_churn["total_rech_amt_8"]))
print("median:")
print(np.median(df_not_churn["total_rech_amt_8"]))
print(np.median(df_churn["total_rech_amt_8"]))
print("Quantiles:")
print(np.percentile(df_not_churn["total_rech_amt_8"],np.arange(0,100,25)))
print(np.percentile(df_churn["total_rech_amt_8"],np.arange(0,100,25)))
In [49]:
#feature importance of total_rech_num_6
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"total_rech_amt_8").add_legend()
Out[49]:
In [50]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_amt_6","total_rech_amt_7").add_legend()
Out[50]:
In [51]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_amt_6","total_rech_amt_8").add_legend()
Out[51]:
In [52]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"total_rech_amt_7","total_rech_amt_8").add_legend()
Out[52]:
1.3) Univariate Analysis on max_rech_amt_6', 'max_rech_amt_7', 'max_rech_amt_8',¶
In [53]:
print("Mean:")
print(np.mean(df_not_churn["max_rech_amt_6"]))
print(np.mean(df_churn["max_rech_amt_6"]))
print("median:")
print(np.median(df_not_churn["max_rech_amt_6"]))
print(np.median(df_churn["max_rech_amt_6"]))
print("Quantiles:")
print(np.percentile(df_not_churn["max_rech_amt_6"],np.arange(0,100,25)))
print(np.percentile(df_churn["max_rech_amt_6"],np.arange(0,100,25)))
In [54]:
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"max_rech_amt_6").add_legend()
Out[54]:
In [55]:
print("Mean:")
print(np.mean(df_not_churn["max_rech_amt_7"]))
print(np.mean(df_churn["max_rech_amt_7"]))
print("median:")
print(np.median(df_not_churn["max_rech_amt_7"]))
print(np.median(df_churn["max_rech_amt_7"]))
print("Quantiles:")
print(np.percentile(df_not_churn["max_rech_amt_7"],np.arange(0,100,25)))
print(np.percentile(df_churn["max_rech_amt_7"],np.arange(0,100,25)))
In [56]:
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"max_rech_amt_7").add_legend()
Out[56]:
In [57]:
print("Mean:")
print(np.mean(df_not_churn["max_rech_amt_8"]))
print(np.mean(df_churn["max_rech_amt_8"]))
print("median:")
print(np.median(df_not_churn["max_rech_amt_8"]))
print(np.median(df_churn["max_rech_amt_8"]))
print("Quantiles:")
print(np.percentile(df_not_churn["max_rech_amt_8"],np.arange(0,100,25)))
print(np.percentile(df_churn["max_rech_amt_8"],np.arange(0,100,25)))
In [58]:
sns.FacetGrid(telecom,hue='churn',size=6).map(sns.distplot,"max_rech_amt_8").add_legend()
Out[58]:
In [59]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"max_rech_amt_6","max_rech_amt_7").add_legend()
Out[59]:
In [60]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"max_rech_amt_7","max_rech_amt_8").add_legend()
Out[60]:
In [61]:
sns.FacetGrid(telecom,hue='churn',size=5).map(plt.scatter,"max_rech_amt_6","max_rech_amt_8").add_legend()
Out[61]:
2) Bivariate Analysis¶
In [62]:
telecom2 = telecom.drop(['mobile_number', 'churn'], axis=1)
In [63]:
# Create correlation matrix
corr_matrix = telecom2.corr().abs()
# Select upper triangle of correlation matrix
upper = corr_matrix.where(np.triu(np.ones(corr_matrix.shape), k=1).astype(np.bool))
# Find index of feature columns with correlation greater than 0.80
high_corr_features = [column for column in upper.columns if any(upper[column] > 0.80)]
print("HIGHLY CORRELATED FEATURES IN DATA SET:{}\n\n{}".format(len(high_corr_features), high_corr_features))
- 36 highly correlated features (> 80%) found in the data set.
- NOTE: We are not dropping these columns manually rather using feature reduction methods to do the job. We worked with RFE, LASSO and PCA.
OBSERVATIONS¶
- Data is higly imbalanced only 5.6% data is churn & remaining is non_churn
- From univariate analysis some relationship is emerging to the churn, but no clear relationship found.
- There are many highly correlated features, found 36 features with more then 0.8 correlation.
Next Steps:¶
- Need to balance the data, apply suitable data-balancing method
- Features are highly correlated, use feture reduction methods before modelling.
Model Building and Evaluation¶
Let's first build a decision tree with default hyperparameters. Then we'll use cross-validation to tune them.
In [64]:
# Importing train-test-split
from sklearn.model_selection import train_test_split
In [65]:
# Putting feature variable to X
X = telecom.drop('churn',axis=1)
# Putting response variable to y
y = telecom['churn']
In [66]:
# Splitting the data into train and test
X_train, X_test, y_train, y_test = train_test_split(X, y,
test_size=0.30,
random_state = 99)
X_train.head()
Out[66]:
In [67]:
# Importing decision tree classifier from sklearn library
from sklearn.tree import DecisionTreeClassifier
# Fitting the decision tree with default hyperparameters, apart from
# max_depth which is 5 so that we can plot and read the tree.
dt_default = DecisionTreeClassifier(max_depth=5)
dt_default.fit(X_train, y_train)
Out[67]:
In [68]:
# Let's check the evaluation metrics of our default model
# Importing classification report and confusion matrix from sklearn metrics
from sklearn.metrics import classification_report, confusion_matrix
# Making predictions
y_pred_default = dt_default.predict(X_test)
# Printing classification report
print(classification_report(y_test, y_pred_default))
In [69]:
# Printing confusion matrix and accuracy
print(confusion_matrix(y_test,y_pred_default))
Plotting the Decision Tree¶
We need the graphviz library to plot a tree.
In [70]:
# Importing required packages for visualization
from IPython.display import Image
from sklearn.externals.six import StringIO
from sklearn.tree import export_graphviz
import pydotplus, graphviz
# Putting features
features = list(telecom.columns[1:])
features
Out[70]:
Python requires library pydot and an external software graphviz to visualize the decision tree.
In [71]:
# plotting tree with max_depth=3
dot_data = StringIO()
export_graphviz(dt_default, out_file=dot_data,
feature_names=features, filled=True,rounded=True)
graph = pydotplus.graph_from_dot_data(dot_data.getvalue())
Image(graph.create_png())
Out[71]:
Hyperparameter Tuning¶
The default tree is quite complex, and we need to simplify it by tuning the hyperparameters.
First, let's understand the parameters in a decision tree.
- criterion (Gini/IG or entropy): It defines the function to measure the quality of a split. Sklearn supports “gini” criteria for Gini Index & “entropy” for Information Gain. By default, it takes the value “gini”.
- splitter: It defines the strategy to choose the split at each node. Supports “best” value to choose the best split & “random” to choose the best random split. By default, it takes “best” value.
- max_features: It defines the no. of features to consider when looking for the best split. We can input integer, float, string & None value.
- If an integer is inputted then it considers that value as max features at each split.
- If float value is taken then it shows the percentage of features at each split.
- If “auto” or “sqrt” is taken then max_features=sqrt(n_features).
- If “log2” is taken then max_features= log2(n_features).
- If None, then max_features=n_features. By default, it takes “None” value.
- max_depth: The max_depth parameter denotes maximum depth of the tree. It can take any integer value or None. If None, then nodes are expanded until all leaves are pure or until all leaves contain less than min_samples_split samples. By default, it takes “None” value.
- min_samples_split: This tells above the minimum no. of samples reqd. to split an internal node. If an integer value is taken then consider min_samples_split as the minimum no. If float, then it shows percentage. By default, it takes “2” value.
- min_samples_leaf: The minimum number of samples required to be at a leaf node. If an integer value is taken then consider - -min_samples_leaf as the minimum no. If float, then it shows percentage. By default, it takes “1” value.
- max_leaf_nodes: It defines the maximum number of possible leaf nodes. If None then it takes an unlimited number of leaf nodes. By default, it takes “None” value.
- min_impurity_split: It defines the threshold for early stopping tree growth. A node will split if its impurity is above the threshold otherwise it is a leaf.
Tuning max_depth¶
Let's first try to find the optimum values for max_depth and understand how the value of max_depth affects the decision tree.
Here, we are creating a dataframe with max_depth in range 1 to 80 and checking the accuracy score corresponding to each max_depth.
To reiterate, a grid search scheme consists of:
- an estimator (classifier such as SVC() or decision tree)
- a parameter space
- a method for searching or sampling candidates (optional)
- a cross-validation scheme, and
- a score function (accuracy, roc_auc etc.)
In [72]:
# GridSearchCV to find optimal max_depth
from sklearn.model_selection import KFold
from sklearn.model_selection import GridSearchCV
# specify number of folds for k-fold CV
n_folds = 5
# parameters to build the model on
parameters = {'max_depth': range(1, 40)}
# instantiate the model
dtree = DecisionTreeClassifier(criterion = "gini",
random_state = 100)
# fit tree on training data
tree = GridSearchCV(dtree, parameters, cv=n_folds,return_train_score=True,scoring="accuracy")
tree.fit(X_train, y_train)
Out[72]:
In [73]:
# scores of GridSearch CV
scores = tree.cv_results_
pd.DataFrame(scores).head()
Out[73]:
In [74]:
# plotting accuracies with max_depth
plt.figure()
plt.plot(scores["param_max_depth"],
scores["mean_train_score"],
label="training accuracy")
plt.plot(scores["param_max_depth"],
scores["mean_test_score"],
label="test accuracy")
plt.xlabel("max_depth")
plt.ylabel("Accuracy")
plt.legend()
plt.show()
We can see that as we increase the value of max_depth, both training and test score increase till about max-depth = 3, after which the test score gradually reduces.
Thus, it is clear that the model is overfitting the training data. Next, let's see how the model behaves with other hyperparameters.
Tuning min_samples_leaf¶
The hyperparameter min_samples_leaf indicates the minimum number of samples required to be at a leaf.
So if the values of min_samples_leaf is less, say 5, then the will be constructed even if a leaf has 5, 6 etc. observations (and is likely to overfit).
Let's see what will be the optimum value for min_samples_leaf.
In [75]:
# GridSearchCV to find optimal max_depth
from sklearn.model_selection import KFold
from sklearn.model_selection import GridSearchCV
# specify number of folds for k-fold CV
n_folds = 5
# parameters to build the model on
parameters = {'min_samples_leaf': range(5, 200, 20)}
# instantiate the model
dtree = DecisionTreeClassifier(criterion = "gini",
random_state = 100)
# fit tree on training data
tree = GridSearchCV(dtree, parameters,
cv=n_folds, return_train_score=True,
scoring="accuracy")
tree.fit(X_train, y_train)
Out[75]:
In [76]:
# scores of GridSearch CV
scores = tree.cv_results_
pd.DataFrame(scores).head()
Out[76]:
In [77]:
# plotting accuracies with min_samples_leaf
plt.figure()
plt.plot(scores["param_min_samples_leaf"],
scores["mean_train_score"],
label="training accuracy")
plt.plot(scores["param_min_samples_leaf"],
scores["mean_test_score"],
label="test accuracy")
plt.xlabel("min_samples_leaf")
plt.ylabel("Accuracy")
plt.legend()
plt.show()
We can see that at low values of min_samples_leaf, the tree gets a bit overfitted. At values > 100, however, the model becomes more stable and the training and test accuracy start to converge.
Tuning min_samples_split¶
The hyperparameter min_samples_split is the minimum no. of samples required to split an internal node. Its default value is 2, which means that even if a node is having 2 samples it can be furthur divided into leaf nodes.
In [78]:
# GridSearchCV to find optimal min_samples_split
from sklearn.model_selection import KFold
from sklearn.model_selection import GridSearchCV
# specify number of folds for k-fold CV
n_folds = 5
# parameters to build the model on
parameters = {'min_samples_split': range(5, 200, 20)}
# instantiate the model
dtree = DecisionTreeClassifier(criterion = "gini",
random_state = 100)
# fit tree on training data
tree = GridSearchCV(dtree, parameters,
cv=n_folds, return_train_score=True,
scoring="accuracy")
tree.fit(X_train, y_train)
Out[78]:
In [79]:
# scores of GridSearch CV
scores = tree.cv_results_
pd.DataFrame(scores).head()
Out[79]:
In [80]:
# plotting accuracies with min_samples_leaf
plt.figure()
plt.plot(scores["param_min_samples_split"],
scores["mean_train_score"],
label="training accuracy")
plt.plot(scores["param_min_samples_split"],
scores["mean_test_score"],
label="test accuracy")
plt.xlabel("min_samples_split")
plt.ylabel("Accuracy")
plt.legend()
plt.show()
This shows that as you increase the min_samples_split, the tree overfits lesser since the model is less complex.
Grid Search to Find Optimal Hyperparameters¶
We can now use GridSearchCV to find multiple optimal hyperparameters together. Note that this time, we'll also specify the criterion (gini/entropy or IG).
In [81]:
# Create the parameter grid
param_grid = {
'max_depth': range(5, 15, 5),
'min_samples_leaf': range(50, 150, 50),
'min_samples_split': range(50, 150, 50),
'criterion': ["entropy", "gini"]
}
n_folds = 5
# Instantiate the grid search model
dtree = DecisionTreeClassifier()
grid_search = GridSearchCV(estimator = dtree, param_grid = param_grid,
cv = n_folds, verbose = 1,return_train_score=True)
# Fit the grid search to the data
grid_search.fit(X_train,y_train)
Out[81]:
In [82]:
# cv results
cv_results = pd.DataFrame(grid_search.cv_results_)
cv_results
Out[82]:
In [83]:
# printing the optimal accuracy score and hyperparameters
print("best accuracy", grid_search.best_score_)
print(grid_search.best_estimator_)
Running the model with best parameters obtained from grid search.¶
In [84]:
# model with optimal hyperparameters
clf_gini = DecisionTreeClassifier(criterion = "gini",
random_state = 100,
max_depth=10,
min_samples_leaf=50,
min_samples_split=50)
clf_gini.fit(X_train, y_train)
Out[84]:
In [85]:
# accuracy score
clf_gini.score(X_test,y_test)
Out[85]:
In [86]:
# plotting the tree
dot_data = StringIO()
export_graphviz(clf_gini, out_file=dot_data,feature_names=features,filled=True,rounded=True)
graph = pydotplus.graph_from_dot_data(dot_data.getvalue())
Image(graph.create_png())
Out[86]:
In [87]:
# tree with max_depth = 3
clf_gini = DecisionTreeClassifier(criterion = "gini",
random_state = 100,
max_depth=3,
min_samples_leaf=50,
min_samples_split=50)
clf_gini.fit(X_train, y_train)
# score
print(clf_gini.score(X_test,y_test))
In [88]:
# plotting tree with max_depth=3
dot_data = StringIO()
export_graphviz(clf_gini, out_file=dot_data,feature_names=features,filled=True,rounded=True)
graph = pydotplus.graph_from_dot_data(dot_data.getvalue())
Image(graph.create_png())
Out[88]:
In [89]:
# classification metrics
from sklearn.metrics import classification_report,confusion_matrix
y_pred = clf_gini.predict(X_test)
print(classification_report(y_test, y_pred))
In [90]:
# confusion matrix
print(confusion_matrix(y_test,y_pred))
Observations and Recommendations¶
Based on the decision tree with the tuned hyperparameters, we can make the following observations and recommendations accordingly. If the customer meets the following criterias, then one is likely to churn
- Case 1
spl_ic_mou_6special calls with less than 0.3 minutes of usagetotal_rech_data_amt_6total data recharge amount less than 0.5 INRlast_day_rch_amt_6number of days since last recharge less than 21.5
Case 2
spl_ic_mou_6special calls with less than 0.3 minutes of usagetotal_rech_data_amt_6total data recharge amount less than 0.5 INRjune_vbc_3gvolume based cost for 3g
Case 3
spl_ic_mou_6special calls with less than 0.3 minutes of usageloc_og_t2t_mou_6local outgoing calls less than 6 minutesmax_rech_amt_6maximum regarche
Case 4
spl_ic_mou_6special calls with less than 0.3 minutes of usageloc_og_t2t_mou_6local outgoing calls less than 6 minutesspl_ic_mou_6special incoming minutes of usage less than 99 minutes
It is recommended that the marketing team should redirect their efforts towards the prevention of churn of these customers.
In [ ]: